'''
Created on Oct 20, 2011

@author: oabalbin
'''

#import os
#import subprocess
import logging
from exome.jobs.base import JOB_SUCCESS, JOB_ERROR, up_to_date,job_done,up_to_date_job
from optparse import OptionParser
from exome.jobs.config import ExomePipelineConfig, ExomeLaneConfig
import exome.gatk_cluster.gatk_realignment_cluster as gra
import exome.gatk_cluster.picard_commands_local as mypicard


def realing_around_indels(jname, configrun, thlane, intervals_to_realign,ncores):
    
    '''This function realign reads around known indels sites'''
    from exome.jobs.job_runner import run_local
    jobrunfunc=run_local
    
    # Parameters:
    extra_mem, num_cores = configrun.gatk_use_mem, configrun.gatk_num_cores
    path_to_gatk, path_to_picard, path_to_sam = configrun.gatk_path, configrun.picard_path, \
                                                configrun.samtools_path
    genomes = configrun.genomes['human']
    ref_genome, snpdb_file, indeldb_file = genomes.gatk_ref_genome, genomes.snpdb, \
                                         genomes.indeldb        
    recal_analysis_outputdir=configrun.gatk_recal_anal_dir
    temp_dir =  configrun.gatk_temp_dir

    
    # Files to create
    indexed_bam_file = thlane.align_bam_file
    realigned_bam_file= thlane.realigned_bam_file
    matefixed_bam = thlane.matefixed_bam_file
    dedup_bam = thlane.markdup_bam_file
    log_dir=thlane.gatk_dir
    
    jobidra=-1
    if up_to_date(realigned_bam_file,realigned_bam_file) and up_to_date_job(realigned_bam_file):
        logging.info("[SKIPPED] Realignment around known indels step. BAM file %s file is up to date" % (realigned_bam_file))
        jobidra=0
    else:
        logging.info("Starting Realignment around known indels BAM file %s" % (realigned_bam_file))
        command = gra.realign_bam_only_knownSites(ref_genome, indexed_bam_file, realigned_bam_file,
                                                       intervals_to_realign, indeldb_file, 
                                                       extra_mem*2, ncores, temp_dir,path_to_gatk)
    
        jobidra = jobrunfunc(jname, command, ncores, cwd=log_dir, stdout='realignIndels.log') #depends.extends(jobidks)            
        if jobidra == 0:
                job_done(dedup_bam)

        
    # Fix Mate information in the Realigned Bam file. It syncronizes the
    # information of mate pairs
    if jobidra != 0:
        raise OSError("command: '%s' returned error code %s" % (' '.join(command), jobidra))
    else:
        if up_to_date(matefixed_bam,matefixed_bam) and up_to_date_job(matefixed_bam):
            logging.info("[SKIPPED] Fix Mates step. Fix Mate BAM file %s file is up to date" % (matefixed_bam))
            jobidfm=0
        else:
            logging.info("Starting Fix Mate step %s" % (realigned_bam_file))
            command = mypicard.fixmateInformation(realigned_bam_file, matefixed_bam, 
                                                                 extra_mem, path_to_picard, temp_dir)
            jobidfm = jobrunfunc(jname, command, ncores, cwd=log_dir, pmem=extra_mem, 
                              deps=jobidra, stdout='fixMate.log')
            if jobidra == 0:
                job_done(matefixed_bam)

    if jobidfm  != 0:
        raise OSError("command: '%s' returned error code %s" % (' '.join(command), jobidfm))
    else:
        # Mark duplicates in the Realigned Bam file. It only mark the duplicates
        if up_to_date(dedup_bam,dedup_bam) and up_to_date_job(dedup_bam):
            logging.info("[SKIPPED] Mark Duplicates step. BAM file with marked duplicates %s is up to date" % (dedup_bam))
            jobidmd=0
        else:
            logging.info("Starting Mark duplicates step %s" % (realigned_bam_file))
            command = mypicard.markDuplicates(matefixed_bam, dedup_bam, extra_mem, path_to_picard)
            jobidmd = jobrunfunc(jname, command, ncores, cwd=log_dir, pmem=extra_mem, 
                              deps=jobidfm, stdout='markDuplicates.log')
    
    # write a file with the done
    if jobidmd == 0:
        job_done(dedup_bam)
        return JOB_SUCCESS
    else:
        JOB_ERROR


def main():
    logging.basicConfig(level=logging.DEBUG,
                        format="%(asctime)s - %(name)s - %(levelname)s - %(message)s")
    
    
    optionparser = OptionParser("usage: %prog [options] ")
    optionparser.add_option("-r", "--config_file", dest="config_file",
                            help="configuration run object")
    '''
    optionparser.add_option("-a", "--analysis_file", dest="analysis_file",
                            help="file with experiment configuration") 
    '''
    optionparser.add_option("-j", "--jname", dest="jname",
                            help="jobname") 
    optionparser.add_option("-l", "--lane_file", dest="lane_file",
                            help="lane object with information about the lane to be processed") 
    optionparser.add_option("-d", "--lane_root_dir", dest="lane_root_dir",
                            help="lane object with information about the lane to be processed") 
    
    optionparser.add_option("-i", "--intervals", dest="intervals",
                            help="intervals around indels around which the realignment is going to be performed") 

    optionparser.add_option("-c", "--ncores", dest="ncores",
                            help="ncores to use in the analysis") 

    (options, args) = optionparser.parse_args()
    
    configrun = ExomePipelineConfig()
    configrun.from_xml(options.config_file)
    thlane  = ExomeLaneConfig()
    #job.from_xml(lane_file, config.output_dir)
    thlane.from_xml(options.lane_file, options.lane_root_dir)

    
    return realing_around_indels(options.jname, configrun, thlane, 
                          options.intervals,options.ncores)

if __name__ == '__main__': 
    import sys
    sys.exit(main())
